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Abstract 

We study the statistical properties of solutions to Burgers' equa- 
tion, vt + vv x = vvxx, for large times, when the initial velocity and its 
potential are stationary Gaussian processes. The initial power spectral 
density at small wave numbers follows a steep power-law Eo(k) ~ \k\ n 
where the exponent n is greater than two. We compare results of 
numerical simulations with dimensional predictions, and with asymp- 
totic analytical theory. The theory predicts self-similarity of statistical 
characteristics of the turbulence, and also leads to a logarithmic cor- 
rection to the law of energy decay in comparison with dimensional 
analysis. We confirm numerically the existence of self-similarity for 
the power spectral density, and the existence of a logarithmic correc- 
tion to the dimensional predictions. 
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1 Introduction 



We study the partial differential equation 

v t + vv x = uv xx (1) 

in the limit of vanishing coefficient v. Here v(x, t) can be taken as the velocity 
of a wave field, a function of one- dimensional spatial coordinate x and time 
t, and v as the viscosity of the medium. 

This equation, first introduced by J.M. Burgers |) as a model of hydro- 
dynamic turbulence, arises in many situations in physics. In the theory of 
wave propagation Burgers' equation ([!]) describes the evolution of nonlinear 
waves in non- dispersive dissipative media [p2 ]. An important example of such 



waves is nonlinear acoustic plane waves [|14[] . When a random force is added 
to the right-hand side it describes surface growth phenomena [10, Burg- 
ers' equation is an effective equation on large scales, when parity-invariance 
holds, and arises as an approximation in description of large-scale structures 



in the Landau- Ginsburg equation NTT]. One further non-trivial application is 



to model the formation and evolution of large-scale structures in the Early 
Universe, within the so-called adhesion approximation []5], |7|, [T7|]. 

From the equation (|l|) the following equation for the potential of the 
velocity, 

v = -fa, (2) 

is derived: ^ 

= -{ipx) 2 + v^) xx . (3) 

By using the Hopf-Cole transformation, ip = 2z/ln#, equation (Q) is trans- 
formed to the linear diffusion equation || 0: 

t = v9 xx . (4) 

The solution of is found by standard convolution with the heat kernel, 
which in the limit v — > may be evaluated by the method of steepest descent. 

The solution of Burgers' equation (|l]) in the limit v —>■ then has the 
form: 

v{x,t) = X ~ a{ - X ' t) (5) 

t 
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where a(x,t) is the coordinate where the maximum of the function 



Gix,a,t)=M")-^-^-- ( 6 ) 
is obtained. The velocity potential is found as 

ip(x, t) = max [G(x, a, t)] . (7) 

In the continuum there are two equivalent ways of computing v(x, t). Either, 
the velocity may be found by taking the derivative of the function ip(x,t) 
from (|7|), using that the velocity is the gradient of the velocity potential 
(H), or equation (H) may be used directly. In a discrete realization the two 
ways are still equivalent, provided we make the correct interpolation between 
the grid points. For long times the difference is not important, and will be 
ignored in the following. 

We use the following definition of the spectral density of energy: 

1 r 00 

E(k,t) = — < v(x + s,t)v(x,t) > e tks ds, (8) 

2ll J-oo 

and we are interested in the case where the initial data is a Gaussian process, 
with continuous spectral support around a dominant wave-number k . We 
want to know how the spectrum changes with time, for sufficiently long 
times. Our approach applies when the initial energy spectrum grows faster 
than k 2 when k is much less than k , and decays quickly when k is larger 
than ko- More to the point, we consider the case where the initial spectrum 
is a power-law at small k 

E {k) ~ \k\ n \k\ « k (9) 

where the exponent n is a real number greater than two. At large k we 
have a smooth cut-off. The construction of the initial conditions including 



normalization factors is described below in section 5.1. 



2 Summary of dimensional predictions 

We assume a dominant initial wave band around ko, which implies a corre- 
lation length of the initial data Zq ~ l/^o- 
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From the characteristic length and the fluctuations of initial velocity we 
can form a characteristic time t n i, the nonlinear time or the characteristic 
time of shock formation: 

t nl = Jg (10) 



<vl > 

As an illustration we show in fig. 1 one typical realization of the initial 
conditions, and how it develops in time. The characteristic scale Iq can be 
estimated as the average distance between local maxima (or minima) of the 
smooth initial data (fig. la). At the nonlinear time the initial data has 
changed to a saw-tooth wave with approximatively the same characteristic 
scale (fig. lb). At larger times the velocity gradient between the shocks 
decreases, and the characteristic scale increases (fig. lc). This behaviour is 
intimately related to the change of the spectral shape with time, to which 
we now turn. 

We will consider in the following only the asymptotic solutions of ([!]) at 
times much greater than the nonlinear time, 

t^t n l. (11) 

The characteristic length scale at this time is derived by balancing in order 
of magnitude the two terms in the Hopf-Cole condition 

2 ^ 72 



— — — < (x - a(x)) > I 

<%> — . (12) 

The distances over which the maximizing is effectively taken at time t are 
therefore such that 

l(t)~(tyj< ^>)\ (13) 
or, in more transparent form, 

l(t)~l Q (t/tra)*. (14) 

We therefore have the picture that at time t the solutions are composed of 
shocks, separated by a typical distance l(t). The spatial gradient of velocity 
between shocks is \ jt from ([5]). The typical velocity and the typical amplitude 
of velocity jumps across shocks is therefore 1(f) /t. From these considerations, 
the mean energy in the velocity field per unit length at time t is: 

E{t) =< v 2 >~ ^ ~< vl > {t/t nl y\ (15) 
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Qualitatively, it is possible to estimate the behavior of the power spectral 
density from the following considerations. As is well known, the development 
of shocks leads to the asymptotic of the spectrum at large wave-numbers, 
E(k,t) ~ k~ 2 , and the parametric pumping of energy to the area at small 
k's leads to the universal quadratic law, E(k,t) ~ k 2 [p], Therefore the 
spectrum should consist of two branches, smoothly transformating from one 
asymptotics k 2 to another k~ 2 . At the turnover point, k*, the spectrum has 
a maximum: 

/ A(t)k 2 , if fc<fc* 
S (*>*)~{ B(t)k- 2 , ifk»k* ■ (16) 



The only possible scale for k* is 



m) ~ ^. (it) 



For the prefactors in (16) , one obtains by comparing with the total energy 



(15) and matching at k*\ 



«S > *S ( ^-V* > (18) 



n2 



B(iH ^(±y 3/2 . (19) 



^0 \tfil 

We see that (jlSf) and (^9j) could also have been derived directly by dimen- 
sional analysis, since A(t) ~ and .B(t) ~ 

It should be noted that if i/ is small but finite we can form a ratio between 
the dissipation time over spatial scale l and the nonlinear time, which we 
can naturally refer to as the initial Reynolds number. 

Reo= it = k<°s> 1 » (20) 

Similarly, at time t we can form the ratio between the dissipation time 
over spatial scale l(t) and the turn-over time, which is on the order of 



l(t)/J< v 2 (t) >, that is, simply proportional to t: 



= w ^ 'q<^> 1/2 (21) 
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Within dimensional considerations the Reynolds number hence remains con- 
stant, because the increase in l(t) cancels with the decrease in typical velocity. 

Another way to introduce the Reynolds number, is to observe that the 
width of a shock is of the order of vj \j< v 2 {t) >. Reynolds number is thus 
also equal to the ratio between the inter-shock distance at time t, that is 
l(t), and the width of the shock. Therefore, as long as Reynolds number 
stays the same, the Hopf-Cole transformation and approximation that the 
solutions have saw-tooth shape are relatively equally accurate at arbitrarily 
large times. 

We will see in the next section that by including logarithmic corrections 
to the dimensional estimates, Reynolds number will in fact slowly decrease 
in time. There will be a characteristic time over which Reynolds number 
decreases to order unity, and then viscosity can not be taken small. This 
time will however be exponentially large in the square of Reo, so if viscosity is 
initially weak it can be taken weak for very long time, although not infinitely 
long. 



3 Beyond the dimensional predictions 

For Burgers' turbulence, in the cases that interest us here, it is possible to 
develop an asymptotic theory going beyond dimensional analysis, which is 
valid at times much greater than t n i. From the mathematical point of view 
we use the property that for large times the parabola in (|5J) varies little over 
a distance of the order of l(t) at its minimum. When l(t) is much larger 
than Zo, and the correlations in the initial conditions decay sufficiently fast 
beyond 1$, then a large number of local maxima of the initial potential ipo(a) 
competes to be the global maximum of the function G(x, a, t), and these local 
minima are practically independent. 

We can then appeal to result in probability theory on the expected max- 



ima of Gaussian processes H, 20]. A more direct, although more cumbersome, 
approach to derive the asymptotics also exists[|J, which has in addition the 
advantage that one may thus also compute order by order corrections to the 
leading behaviour 0. 

We begin by observing that the mean square velocity and mean square 
velocity potential are both finite for the process we consider. By comparing 
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them we can thus form a characteristic length 



<^ 2> (22) 

which will of course also be of the order 1/ko, but with a well-defined finite 
prefactor depending on the functional form of the spectrum. 

Let us consider the decay of energy with time. Averaging both sides (H) 
over space leads to 

d t < -0 >= E{t) (23) 

The decay of energy can thus be derived from the expectation value of the 
velocity potential at time t. The maximization position (a) in ([?]) can at 
most be a distance of leading order \/i away from the Eulerian coordinate 
(x). The expectation value of the velocity potential at time t should therefore 
be similar to the expected maximum value of the initial velocity potential 
on a stretch of the line of length y/t. One would expect this value to grow 
in time, although not very fast. In fact, a well-known result on stationary 



Gaussian processes [H, EO, B] states that this maximum has a doubly exponen- 



tial distribution, such that the expected maximum on a stretch of the line of 
length L is < i/jq > x I 2 J2 log , where < ipQ > has the same meaning as 



introduced above, and l as in equation ([22]). 
Substituting lo(t/t n iY^ 2 for L we arrive at 



log ^- (24) 



Although we should not a priori expect that (|4]) is correct up to constants, 
it turns out to be the case JT2], f|, §, [£| . By differentiating in t follows Kida's 
law for the the deacy of energy at long times: 

E(t) = \<vl> {t/traY 1 (log(t/2nt nl ))- 1/2 (25) 
Again, this expression holds in the asymptotic theory up to additive correc- 



tions that become small relative to (|25|) in the limit 



We can also estimate the behaviour of the characteristic length. When 
time grows, the expected maximum of the term ipo in (0) will slowly become 
larger than (x — a) 2 /2t, if the effective length of the maximization operation, 
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(x — a), is just taken to be of order The balance equation, (|I2"D, therefore 
slowly becomes inappropriate to fix the spatial scale. Let us therefore instead 
assume that the maximum is attained at a value of (x — a) around L. The 
attained value of if)(x,t) would then typically be 

< ^ >l ~ (V< i>l >j21og— O - ^. (26) 



2t 



If we now maximize ( |26| ) over L we find that the maximum is attained at a 
value of 

/(t)~/o(t/U 1/2 Qog^-)^ (27) 

27rt n ; 

The characteristic length-scale therefore grows slightly more slowly than the 
dimensional estimate. Let us note that the estimate of the maximum (f2~3|) 
remains unchanged. Only subleading terms are influenced by the small cor- 
rection to the characteristic length. Let us also note that the dimensional 
estimate E(t) ~ l 2 (t)/t 2 remains correct up to logarithmic corrections, which 
indicates that the fluctuations are very small. 

In a similar spirit one may also compute two-point probability densities 
[|, ||, correlation functions and the power spectral density ]12|, |, ||, ||. In 
particular, the energy spectrum has the expression 

l 3 (t) - 

E(k,t) = ^E(kl(t)), (28) 



where l(t) is the characteristic scale in (p7|) , and the function E(k) can be 
computed to the same leading order in the asymptotics^, [| |12| . 

In fig. 2 (a) three power spectra (averaged over 1000 realizations of the 
random process) are shown at three different moments of time. The initial 
spectrum was k 3 at small k. The preservation of the shape of each curve 
is evident. It is also seen that the total energy (the area under the curves) 
decreases, and so does the characteristic wave- number k*. Fig. 2 (b) contains 
three power spectra computed at the same time at the self-similar stage 
(averaged over 1000 realizations). The initial spectra were k 3 , k 6 and k 12 at 
small k, respectively. 

The analytical expression for E is somewhat unwieldy. However, qualita- 
tively, E is simply described as a smooth interpolation between two power- 
law asymptotes that hold when k is respectively much smaller or much larger 
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than unity. In the first case E(n) ~ 1.08/t 2 , in agreement with the expected 
leading k 2 behaviour of the spectrum according to (|T6|) , while for large k 
E(k) ~ 0.359/k 2 , also in agreement with (JT^) []5J. In Fig.3 E(k) is shown 
together with numerical results obtained by averaging over 1000 realizations 
of the initial conditions. The numerically obtained curve (as a result of av- 
eraging over 1000 realizations of the initial random process) is artificially 
shifted slightly downwards to be distinguishable from the analytical curve. 

It is seen that the asymptotic formula describes the numerical data very 
well, not only in the limits of relatively large and small wave numbers, but 
also at the top, where the spectrum switches between the two asymptotes. 
From fig. 3 is also seen that the transformation region from the top down 
towards k~ 2 is thin, but the transformation region from the top down to k 2 
is considerably wider. 

In fig. 4 the function l(t) is shown. The dashed line corresponds to the 
calculation based on the expression fll3"|) . The result of the computations 
taking into account the logarithmic correction (|27|) is shown by the solid line. 
The centered symbols show results from numerical simulations of an ensemble 
of 500 realizations. The length scale was estimated from the numerical data 
in the following way: the averaged spectrum was first found, and by fitting 
the two asymptotes Ak 2 and Bk~ 2 , their intersection was taken as the main 
wave number k* . The numerical estimate of l(t) was then defined as 1/k*. 
The fact that the results with the logarithmic corrections fit the data better 
proves numerically the validity of the asymptotic theory. 

In fig. 5 we show the results of a procedure analogous to that depicted in 
fig. 4 but for the function E(t). It is seen again that the solid line obtained 
by taking into account the logarithmic correction fits the numerics (centered 
symbols) better. Since the correction is slow, the difference between this 
result and the result of the dimensional analysis (dashed line) is visible but 
small. 

From which have been shown above one might conclude that the logarith- 
mic correction is not very significant, since it gives a small addition to the 
dimensional predictions. To a certain extent that depends on what quantity 
one looks at. We can for instance compute the prefactors of the power-law 
for the spectrum at large and small wave numbers. If E(k) ~ A(t)k 2 up to 
a maximum wave number k*, then the total energy at time t is to leading 
order A(t)(k*) 3 . But we know already the energy and k* as functions of t, so 
by comparing we compute A(t). By matching at k* we then also find B(t), 
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and the results are: 



w**>s(B 1 ' (29) 

In the case of A(t), the correction to the dimensional predictions is more 
significant than those shown previously. In fig. 6 the centered symbols show 
the values of A(t) obtained from the numerical simulation. The dashed line 
is the dimensional prediction computed according to (pT8|), while the solid line 
takes into account the logarithmic correction (|29|). 

There is also a real physical interest in the logarithmic corrections, since 
they change qualitatively the behaviour at large times. Let us form the 
Reynolds number at time t, as in (f21|): 



Re(t) = WV ~ Re (log-?-) (31) 

where Re® is the initial Reynolds number as in (p0|) . Although we compute 
the solutions to Burgers' equation in the limit where viscocity vanishes, that 
is, in the infinite Reynolds number limit, Reynolds number still decreases 
with time [|1], [5], |!J. Sooner or later we will reach a stage where viscosity 
can not be taken arbitrarily small, and at even later time we will eventually 
reach a the linear stage of decay. The time to reach the linear regime can be 
estimated by taking Re(t[) ~ 1, which leads to 

ti ~27rt nl exp{Re 2 ). (32) 



4 Different dimensional ansatze 

We started our more precise analysis by expressing the correlation length 
as Iq = J • We might also consider the initial velocity gradient, and 
its variance, that we write o\ =< (d x Vo) 2 >. We could then define the 

correlation length by l' Q = \j <v ^ , and a new nonlinear time by t' nl = l/a u . 
This definition of the nonlinear time has some appeal, since a shock first 
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forms where the velocity gradient has a local minimum, and does so after a 
time l/|(9a;fo|. 

For initial spectra going as k n with n greater than two, 1' and t' nl are of 
course simply proportional to 1$ and t n [. In dimensional analysis it does not 
matter which ones we use, but, as we have seen, in the asymptotic analysis the 
correlation length and the nonlinear time that appear are those given above 
and used in section |^. Furthermore, a considerable part of the asymptotic 
analysis is valid down to n equal to one [§]. In that limit Iq diverges, while 
Vq remains finite. The expressions from dimensional analysis using 1' and t' nl 
therefore grow progressively less accurate as n tends to one from above. 

If we lower n further and consider the interval [—1, 1], Iq and t n \ do not 
exist, because the initial velocity potential is then not a stationary process, 
only a process with stationary increments. The initial velocity will still be 
stationary however, and its correlation length will be approximatively l' Q , and 
the shock formation time will still be about t' nl . For distances \a — x\ much 
greater than 1' , < (ipo(a) — ipo(x)) 2 > is equal to (3 2 {\a — x|) 1_n , where f3 2 is 

2 

given as a 2 r ( 2 _ TO 4 ) 7r cos ™ • Balancing in (|7|) gives l(t) ~ (/3t)z+^. We see that 
(3 2 diverges as n tends to one from below. 

If we instead proceed by dimensional analysis using 1' , then < (ipo(a) _ 
ip{i{x)) 2 > should be similar to (1' ) 2 < v I > (^jr^) 1 ^ n , with a numerical 
factor which depends on how the cut-off is made at large k. Balancing 

2 

in (|7|) gives l(t) ~ l'o(t/t'ni) 3+n ■ The upshot of this discussion is that t' nl 
is a good measure of shock formation time, and can be used as a nonlinear 
time in the dimensional estimates in the whole interval of n greater than — 1, 
except around one, where numerical prefactors diverge, because the long- 
term asymptotics changes qualitatively around this point. 

5 Numerical work 

5.1 Normalizations of initial spectrum 

We use the following smooth cut-off of the initial power spectrum: 

_ k 2 

E {k) = a 2 n \k\ n e (33) 
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We let k vary in the first Brillouin zone, that is k goes between —1/2 and 
1/2, with ko of the order of 1/2. The initial spatial scale l is therefore of 
order one. 

The variances of the velocity and the velocity potential can be computed 
to be: 

„2 ^ _ ^,2 F ( , n +l 



<«S> = k5+\ (34) 

<$> = a£r(=^i)*r\ (35) 

where T(x) is the gamma-function. 

Furthermore, we use the convention that < ip^ > is equal to one, which 
is of course equivalent to making a choice of scale of time. We have then 

<<> = (36) 



and, following the above, 



^0 = 7- 



r 



n-l 
2 



and 



A:o^r(H±i) 

r(s=i) 



(37) 



(38) 



nl r(=±i)*g- 

We can also compute the variance of the velocity gradient, which is 

(39) 

We can hence form the second estimate of the initial correlation lengths, 
Zq, and the ratio between the two lengths, which is 

l'o _ r (¥) 



(40) 



We see indeed that these two lengths only differ by a constant depending on 
n, as was stated in the previous discussion. 
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5.2 Generation of initial conditions 



Fourier components of a Gaussian process are independent Gaussian vari- 
ables. We therefore synthesize the initial potential of the velocity by first 
generating random variables distributed according to 



P(a k ) = 




where 

a\ = E^ {k)dk. 

Here is used the well known relation between the power spectra of the process 
and its derivative 

E^ {k) = k- 2 E {k) (41) 

where the form of Eo(k) is chosen with a smooth cut-off at large k according 
to (|33|) . By an inverse Fourier transform we find the initial potential in real 
space, and repeating the whole process many times we sample the desired 
ensemble of Gaussian initial conditions. 



5.3 Fast Legendre Transforms 

In numerical simulations the initial data are always generated as a discrete 
set of N points. It could be assumed naively that the number of operations 
necessary to compute the maximization ([?]) for all values of x scales as 0(N 2 ). 
It may however be shown, using (|7]) that a(x) is a nondecreasing function of 
x. The number of operations needed in an ordered search therefore scales 
as 0(iV log 2 N) when using the so-called Fast Legendre Transform procedure 

MM- ' 

When the number of grid points is large they can not all fit into the 
working memory of the CPU at the same time. Even if the operations count 
is the same, it is preferable to limit paging of data from an external memory 
storage to a minimum. This is achieved in the in-order algorithm of Noullez 



[fL3|l that we use here. 



5.4 Discrete Fourier analysis 

After Legendre transform, we have at our disposal a function ip{x, t) sampled 
at iV discrete points ip(xi,t), 1 — 1,2, ...,N. The sample interval is Ax = 1. 
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Computing the discrete Fourier transform by using FFTs as 

^(j,t)=^(xi,t)e 2 ^ N j = 0,...,N-l (42) 
i 

we will get Fourier coefficients at wave numbers kj = ^F, where the last half 
indicates negative frequencies. Since the signal is real the Fourier amplitudes 
are related by if)(kj,t) = if)*(—kj,t). 

Knowing the spectral density of the potential ip(kj,t), it is easy to find 
the spectral density of the velocity, since 

v(kj,t) = ikjip(kj,t). 

Different definition can then be made of the numerically obtained power 
spectrum, which only differ by prefactors of N. We use here the normaliza- 
tions given in a standard reference |i~5| : 



E(o,t) = 4^(MI 2 



N 2 

E(k h t) = A|«(fc,,t)| a , j = l,2,...,(y-l), (43) 
E{k c ,t) = -L\v(k c ,t)\ 2 . 

Since the initial potential ipo(xi) is obtained by an inverse Fourier trans- 
form it is a periodic function with period N. But the use of the fast Legendre 
transform of a periodic function introduces undesirable edge effects which do 
not allow to achieve periodicity in the resulted potential i/j(xi, t). This causes 
distortion of the spectrum connected with the finite, but not periodic, func- 
tion in real space. We eliminate this effect by constructing the periodic 
potential at time t for the cost of using the fast Legendre transform over ad- 
ditional N points of the initial potential. The 'new' periodic initial potential 
in obtained by simple addition a fragment of ipo(xi) with I in [N/2 + 1,N] 
to the left-hand side of the same realization. A fragment with I in [1, N/2] 
is then added to the right-hand part. Doing this we again obtain a periodic 
function ip' (xi) but now with 'double' periodicity 

ip'o(xi) = 4>o(xi +N/2 ), l<l<N/2, 

^/ (xi) = Mxi-N/2), N/2 + l<l<N + N/2, 

Vq(xi)= M*i-n-n/2), N + N/2 + l<l<2N. 
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After performing the fast Legendre transform over 2N points we get a func- 
tion ip'(xi,t). This function loses the 'double' periodicity, but the middle N 
points of this realization 

^(x h t)=^'(xi +N/2 ),l< I < N 

form a periodic function ip(xi, t) with period iV which we use further for the 
Fourier analysis. 

We find that to check the described theory not very large simulations are 
needed. In Fig.l - 6, shown above, we have used N equal to 2 15 , i.e. 32768. 
No qualitative difference appears with the use of larger N or greater number 
of realizations. 

6 Acknowledgments 

This work was supported by the Swedish Natural Science Research Coun- 
cil through grant S-FO-1778-302 (E.A.), by the KTH CITEC Project No. 
930-131-95 (S.N.G.), and by the Swedish Institute (S.I.S.). E.Aurell and 
S.I.Simdyankin thank the Center for Parallel Computers for hospitality. 

References 

[1] Burgers J.M., 1974, "The Nonlinear Diffusion Equation", D. Reidel, 
Publ. Co.. 

[2] Cole J.D., 1951, "On a quasi-linear parabolic equation occuring in aero- 
dynamics", Quart. Appl. Math., 9, 225. 

[3] Cramer H. and Leadbetter M.R., 1967, "Stationary and Related 
Stochastic Processes" , New York: Wiley. 

[4] Fournier J.D. and Frisch U., 1983, "L'equation de Burgers deterministe 
et statistique" , J. de Mec. Theor. et Appl, 2, 699. 

[5] Gurbatov S.N., Malakhov A.N., and Saichev A.I., 1991, "Nonlinear ran- 
dom waves and turbulence in nondispersive media: waves, rays, parti- 
cles", Manchester University Press. 



15 



[6] Gurbatov S.N. and Saichev A.I., 1981, "Degeneracy of one-dimensional 
acoustic turbulence at large Reynolds numbers", Sov. Phys. JETP, 53, 
347. 

[7] Gurbatov S.N., Saichev A.I. and Shandarin S.F., 1989, "Large Scale 
Structure of the Universe in the frame of the model equation of nonlinear 
diffusion", Mon.Not.R.Astr.Soc, 236, 385-405. 

[8] Gurbatov S., Aurell E., Frisch U., Simdyankin S. and Toth G. (in prepa- 
ration) . 

[9] Hopf E., 1950, "The partial differential equation u t +uu x = u xx " , Comm. 
Pure Appl. Mech., 3, 201. 

[10] Kardar M., Parisi G. and Zhang Y.C., 1986, "Dynamical scaling of grow- 
ing interfaces", Phys. Rev. Lett., 56, 889. 

[11] Kuramoto Y., 1985, "Chemical Oscillations, Waves and Turbulence", 
Springer Verlag, Berlin. 

[12] Kida S., 1979, "Asymptotic properties of Burgers turbulence", J. Fluid 
Mech., 93, 337. 

[13] Noullez A. and Vergassola M., 1994, "A fast algorithm for discrete Leg- 
endre transforms" J. Sci. Comp., 9, 259-281. 

[14] Rudenko O.V., Soluyan S.I., 1977, "Theoretical foundations of nonlinear 
acoustics" , New York: Plenum Press. 

[15] Press W., Flannery B., Teukolsky S and Vetterling W. 1988, "Numerical 
Recipies in C", Cambridge: Cambridge University Press. 

[16] Schwartz M., Edwards S.F., 1992, "Nonlinear deposition: a new ap- 
proach", Europhys. Lett., 20, 301-305. 

[17] Shandarin S.F. and Zel'dovich Ya.B., 1989, "The large-scale structure of 
the Universe: Turbulence, intermittency, structures in a self-gravitating 
medium", Rev. Mod. Phys., 61, 185. 

[18] She Z.S., Aurell E. and Frisch U., 1992, "The inviscid Burgers equation 
with initial data of Brownian type" Comm. Math. Phys., 148, 623-641. 



16 



[19] Sinai Ya.G., 1991, "Statistics of shocks in solutions of the inviscid Burg- 
ers equation", Comm. Math. Phys., 148, 601. 

[20] Tikhonov V.I., 1970, "Vybrosy sluchainykh protsessov (Exceedances of 
Random Processes)", Moscow: Nauka 

[21] Vergassola M., Dubrulle B., Frisch U. and Noullez A., 1994, "Burg- 
ers' equation, devil's staircase and the mass distribution for large scale 
structures" Astro. Astrophys., 289, 325-356. 

[22] WhithamG.B., 1974, "Linear and Nonlinear Waves" , New York: Wiley. 



17 



Figure captions: 



Figure 1: A realization of the velocity field at time t — (a), evolved process at 
dimensionless time tjt n i ~ 10 5 (b), and the same realization at a larger 
time t/t n i = iO 6 (c). 

Figure 2: (a) Power spectral densities computed at different times: (1) t/t nl ~ 
5 • iO 8 , (2) t/t nl ~ 3 • iO 9 , (3) t/t nl ~ 2 • 10 10 . Log-log plot, 
(b) Three power spectral densities computed at the same time. The 
exponent n of the initial power spectrum for the curve above n — 3, 
and for the other two, which are hardly distinguishable, n = 6 and 
n = 12. Log-log plots. 

Figure 3: Power spectral densities resulted from evaluation of an analytical inte- 
gral expression (upper curve) and from the direct numerical simulation 
(lower curve shifted down to be distinguishable from the upper one). 
Log- log plot. 

Figure 4: Characteristic space-scale as a function of time. The prediction l(t) ~ 
\ft (dashed line) is by dimensional analysis, while the more precise 
approximation is (logt) -1 / 4 \/t (solid line). Log-log plot. 

Figure 5: Total energy as a function of time E(t). 

Centered symbols - the numerically obtained values of E(t). The pre- 
diction E(t) ~ 1/t (dashed line) is by dimensional analysis, while the 
more precise approximation is 1 / ty/\og t (solid line). Log-log plots. 

Figure 6: The asymptotic prefactor A(t) at the IR part of the power spectral 
density as a function of time. Centered symbols - values of A(t) re- 
sulted from the numerical simulation. Dashed line - the dimensional 
prediction, and solid line - analytical result with taking into account 
the logarithmic correction. 



18 



